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We present a novel method for estimating the surface horizontal velocity on ice shelves using laser altimetry 
data from the Ice Cloud and land Elevation Satellite (ICESat; 2003-2009). The method matches undulations 
measured at crossover points between successive campaigns. Elevation measurements are first relocated 
into a time-varying (moving) coordinate system using an initial velocity (e.g., from VELMAP), and then cross- 
over height differences are minimized with an adjustment vector. Errors in geolocation of the ICESat tracks 
result in some error in the adjustment vectors, but these are small relative to the velocity adjustment for 
fast-moving ice shelves. We use the algorithm to estimate changes in the ice velocity of Ross Ice Shelf be- 
tween an earlier mapping (from VELMAP) and the ICESat period. The new velocity field is compared with ve- 
locities from in situ measurements and satellite radar interferometry. The slowdown of 98 ± 34 m yr“ 1 
(-23%) is observed in the ice shelf downstream of Whillans Ice Stream, and the deceleration rate is 3.1 ± 
1.1 m yr~ 2 during last three decades. The method can be expanded to the simultaneous mapping of ice hor- 
izontal velocity, ice thickness change, and surface deformation for Antarctic ice shelves as well as a more 
accurate mapping using future ICESat-2 measurements. 

© 2012 Elsevier Inc. All rights reserved. 


1. Introduction 

The ice velocity of an ice shelf is a key parameter for understanding 
its dynamics and its interaction with the inland ice sheet. A changing ice 
shelf velocity field provides insight into dynamical changes of the 
inflowing grounded ice, or changes in basal melt or stress balance with- 
in the shelf. Many methods exist for mapping ice shelf velocity (e.g., in- 
terferometric SAR (InSAR) or InSAR speckle tracking; also, image 
feature tracking using satellite visible-infrared images) (Bindschadler 
et al., 1994; Joughin, 2002; Rignot et al., 2008; Scambos et al., 1992). 
The application of InSAR is most accurate but a tracking of long-term 
ice velocity is not available in fast-moving areas because of 
unrecoverable unwrapped phases. The intensity correlation of repeat 
visible-infrared satellite images (e.g. Landsat images) is one of the com- 
mon methods to measure surface velocity albeit with large time steps. 
However, it relies on optical surface contrast features that are naturally 
scarce over ice shelves, and a lack of ground control point on large ice 
shelves causes large geolocation errors. Both methods are in general 
limited by data availability, particularly at the highest latitudes of the 
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ice shelves. In contrast, ICESat provides improving geolocation accuracy 
and coverage in these areas, with increasing crossover density up to the 
satellite latitude limit of 86°S. 

We present a new method to estimate the ice velocity of the large 
Antarctic ice shelves using the elevation measurements from the Ice 
Cloud and land Elevation Satellite (ICESat). The Geoscience Laser Altim- 
etry System (GLAS) onboard the ICESat measures elevation with a 
footprint size of -65 m at a space interval of 172 m. Elevation measure- 
ment precision is about a decimeter, although under ideal conditions a 
few centimeters can be achieved (Fricker et al., 2005; Schutz et al., 
2005). The accuracy and high resolution of ICESat measurements make 
it possible to use satellite altimetry as a method of estimating surface 
velocity on ice shelves. 

Abdalati and Krabill (1999) have previously determined ice veloc- 
ities in portions of the Jakobshavn drainage basin using airborne laser 
altimetry data. Their technique is based on dense measurements 
enough to be interpolated onto regularly spaced grids suitable for 
cross-correlation analysis but densely scanned laser altimetry data 
are rarely available for large ice shelves. More recently, the correla- 
tion between two repeat ICESat tracks have been used to estimate 
the ice shelf surface velocity (Marsh & Rack, 2012). However, this 
method cannot be used as a standalone method as it currently 
measures only the along-track surface velocity. Instead of these 



252 


C-K. Lee et al. / Remote Sensing of Environment 124 (2012) 251-258 


correlation analyses, we refine the ice velocity measurement at each 
crossover location by minimizing the elevation difference of the 
crossovers re-located in the moving coordinate frame (movement 
prescribed by the velocity field). Below, we describe this algorithm 
in detail and then present a revised ice velocity of the Ross Ice Shelf 
(RIS) for 2003-2009 periods. 

2. Data 

2.1. ICESat data 

ICESat data products GLA12 (release 531, for ice sheet elevations) op- 
erated from October 2003 to April 2009 (15 campaigns) are used in this 
study. We apply an inverse barometer correction to remove the atmo- 
spheric loading using the European Centre for Medium-Range Weather 
Forecasts' ERA-Interim data (Simmons et al., 2007) along with the satu- 
ration correction. Our analysis of various tide models over the RIS finds 
that the ICESat measurements corrected by the Circum-Antarctic Tidal 
Simulation model (CATS2008a), an update of Antarctic regional models 
by Padman et al. (2002) yields an RMS crossover height difference of 
-0.17 m, whereas those corrected by the GOT99.2 model, which is al- 
ready applied to current ICESat data products, yields -0.27 m. Therefore, 
we use the CATS2008a to substitute for GOT99.2. 

2.2. Crossovers 

Crossovers are the intersections of ascending and descending 
ICESat laser profile tracks. We focus particularly on the inter- 
campaign crossovers (crossovers between tracks from two different 
campaigns) so that a sufficient time has passed between the surface 
elevation measurements to refine the velocity. Since we use all ICESat 
profiles acquired in 15 campaigns, a maximum of 210 inter-campaign 
crossovers are available for a crossover point. The intra-campaign 
crossovers (crossovers between tracks from a campaign) are used to 
assess the measurement errors in ICESat elevation profiles because 
the effects from ice flow are negligible within a campaign (the 
mean time difference of intra-campaign crossovers is- 11 days). The 
crossover measurement density significantly increases toward the 
ICESat maximum latitude (~86°S) (Fig. 1). For example, the number 
of inter-campaign crossovers from 15 ICESat campaigns in an area 
of 0.5° (latitude) x 2° (longitude) is -500 at the latitude of 78°S, 
while it exceeds 2000 at the latitude of 84°S. 



Fig. 1 . ICESat tracks of Laser 2a (yellow lines) and velocity observations from the 
VELMAP (blue circles) used to set the initial velocity in the RIS. The red squares are 
the GPS stations occupied in 2005 (Brunt, 2008), and the red dots are three RIGGS sta- 
tions for validation mentioned in Section 4.3. 


2.3. Velocity data for the Ross Ice Shelf 

The Antarctic Ice Velocity Data (VELMAP) distributed by National 
Snow and Ice Data Center (NSIDC; http://www.nsidc.org/data/velmap) 
are used as the initial velocity fields. The VELMAP data set in the RIS is 
a compilation of several sources acquired at several times from 1970s 
to the 1990s (Fig 1). The Ross Ice Shelf Geophysical and Glaciological 
Survey (RIGGS) obtained the ice velocity by the point-positioning meth- 
od from 1972 to 1974 over the entire RIS (Thomas et al., 1984). Along 
the Siple Coast and downstream of Skelton Glacier, many velocity data 
were acquired by the correlation method of Landsat images from 1987 
to 1990 (Scambos et al., 1992). The velocities of Byrd Glacier were deter- 
mined by photogrammetric methods from two sets of aerial photo- 
graphs taken in 1978 and 1979 (Brecher, 1986). 

3. Method 

3.1. Relocation of ICESat profile into the moving coordinate frame 

Two sections of laser elevation profile at an inter-campaign cross- 
over do not refer to the same location on the ice surface due to flow of 
the ice sheet. In order to find the crossover area on the moving ice, we 
relocate the ICESat measurements into a new coordinate system mov- 
ing with a velocity v. If the origins of the stationary (geodetic) coordi- 
nate system and the moving coordinate system coincide at time t 0 
(reference time), the position vector in the moving coordinate sys- 
tems (r) is described under the Galilean transform: r = r— v(t— 1 0 ), 
where r is the position vector in the geodetic coordinate system and 
t is the measurement time. We then find new sets of crossovers dis- 
placed from the crossovers in the geodetic coordinate system. Hereaf- 
ter, we call the conventional crossover in the geodetic coordinate 
system “orbital crossover” and the new crossover in the moving coor- 
dinate system “ ice-surface crossover” (Fig. 2). The elevation differ- 
ence between two shot points at the ‘ice-surface’ crossover will 
approach zero if velocity of the moving coordinate frame is valid 
and elevation changes from other effects (i.e. snow accumulation, 
basal melting, and spreading) are negligible (Fig. 2). 

3.2. Factors influencing the crossover differences 

The crossover elevation differences of orbital crossovers (Ah) con- 
tain elevation differences due to ice flow (Ahf), elevation difference 
due to ice thickness change (Ah t ), and measurement error (e m ): 
Ah = Ah { + Ah t + e m . The elevation differences caused by ice flow 
(Ahf) play the key role in our derivation of the velocity fields. Consider- 
ing the mean ice velocity in the RIS (-500 m yr _ ^ and the mean time 
difference of inter-campaign crossovers (-2 yr), the amount of ice 
movement inherent in inter-campaign crossovers is -1 km in average. 
The variability of Ah { is thus expected to be about 0.5 m roughly mea- 
sured by the standard deviation of elevation differences at 1 km lag. On 
the other hand, the ice thickness change (Ah t ) is mainly dependent on 
the snow accumulation (Ahf), basal melting (Ahf), and deformation 
(Ahf). Compilations using in situ data and the Advanced Microwave 
Scanning Radiometer (AMSR-E) show that the accumulation rate 
ranges between 106 to 147 kg m -2 yr -1 in the RIS, which was 
converted to a firn-equivalent accumulation rate (Ahf of 0.30 to 
0.42 m yr - 1 using a firn density of 350 kg m -3 (Arthern et al., 2006). 
Glacial melt from the RIS was estimated to be 33-55 km 3 yr - 1 from 
the data-based analysis of melt water production (Loose et al., 2009). 
Applying the area of RIS (-487,000 km 2 ), the basal melt rate (Ahf) in 
the overall RIS is about 0.07 to 0.11 m yr -1 . If we assume that the 
strain rate in the RIS is in the range from 0.0005 to 0.001 yr -1 (Alley 
et al., 2008), the thickness change rate due to ice deformation (Ahf) 
is approximately —0.8 to — 0.15 m yr -1 using the strain rate correc- 
tion, which is defined as the product of strain rate and ice thickness. 
Combining all three effects (Ah a — Ah\> + Ahf), the total thickness 
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Fig. 2. Comparison between the orbital and ice-surface crossovers of track 44 and 1288 at 78.5°S and 177.7°W. Left panel: the original ICESat tracks in the geodetic coordinate sys- 
tem (blue lines) are relocated in the moving coordinate system with the reference time of January 1, 2004 using a velocity v (red lines). The yellow circles indicate the ice-surface 
crossovers along the Laser 2a ascending track. The background image is the MODIS mosaic of Antarctica (MOA, Nov. 2003-Feb. 2004). Right panels: the elevations along a Laser 2a 
track (black lines) are compared with the elevations from other tracks (yellow circles) at the ice-surface crossovers according to the initial (v 0 ) and adjusted (v) velocity. The blue 
circles indicate the elevations at the orbital crossovers. 


change rate varies in the range from — 0.6 to 0.2 m yr _ \ Since the el- 
evation change is about 0.1 1 times the thickness change of ice in sea 
water, the Ah t during the mean time difference of inter-campaign 
crossovers (-2 yr) is in the range from —0.14 to 0.04 m taking up a 
portion much smaller than the effect of ice flow (~ 0.5 m). 

The measurement error (e m ), which includes the errors in the sen- 
sor, atmospheric and tide correction, wind transport of surface snow, 
crossover interpolation, etc., can be assessed by the elevation difference 
of the intra-campaign crossover (crossover between two tracks from a 
campaign) because the effects from ice flow and ice thickness change 
are negligible within a campaign (the mean time difference of intra- 
campaign crossovers is -11 days). The RMS elevation difference of 
intra-campaign crossovers is -0.17 m over the RIS. The geolocation 
error of the ICESat tracks result in some error in the velocity estimation, 
but these are small relative to the ice movement (-1 arcsec, or -3.5 m) 
(Luthcke et al., 2005). Consequently, since the ice flow is the most im- 
portant factor influencing the crossover elevation difference and the 
measurement error is random, we expect a reasonable estimate of the 
ice velocity by means of searching on the velocity of the moving coordi- 
nate frame, seeking to minimize the elevation differences of ice-surface 
crossovers (to a limit approaching -0.17 m on average over the shelf). 

3.3. Velocity refinement by minimizing the crossover difference 

We apply the RMS elevation difference of ice-surface crossovers 
to the fitting function (f) in the minimization formula, expressed as 

minf = min ^E [h' a ( r'j) — h' d (r' f )] 2 (1 ) 

where h a (r'^ and h d (r'^j are the elevations measured along ascend- 
ing and descending profiles at the position of the ith ice-surface 
crossover (r' f ) respectively. 

We first examined the behavior of our fitting function in a test area 
(Fig. 3). The fitting functions were exhaustively calculated in a pre- 


defined velocity space with a maximum deviation of 200 m yr - 1 from 
a given initial velocity and an interval of 1 0 m yr - 1 (Fig. 3b). The veloc- 
ity deviation from the initial velocity can be determined by finding the 
position of minimal fitting function (100 m yr -1 in an eastward direc- 
tion and 10 m yr -1 in a northward direction in this example). Fig. 3c 
and Fig. 3d compare the elevations along the ICESat tracks relocated 
into the moving coordinate systems with the initial velocity and the ve- 
locity estimated from the minimal fitting function, respectively, and 
superimposed on a MODIS composite mosaic of Antarctica (MOA, 
Nov. 2003-Feb. 2004) (Scambos et al., 2007). Looking into three black 
circles, the features in the ICESat tracks relocated with the estimated ve- 
locity are remarkably well aligned with the linear features in the MODIS 
mosaic, indicating the improvement of ice velocity. 

The quadratic shape and unique minimum of fitting function indi- 
cate that the solution can be efficiently found by the non-linear least 
square inversion. We use the trust-region non-linear least square 
method with two unknown parameters (horizontal velocity v x and 
v y ) in order to minimize the fitting function. The trust-region method 
is a kind of iterative methods for optimization, which is robust and 
can be applied to ill-conditioned problems, similar to the classical 
Levenberg-Marquardt method. In the trust-region algorithm for 
non-linear least squares, the approximate model such as quadratic 
approximation is only trusted in a region near the current iterate so 
that the update of solution is bounded in the “trust region” updated 
from iteration to iteration (More & Sorensen, 1983). 

4. Result 

4.1. Velocity adjustment over the RIS 

To increase the signal to noise ratio, the area patch, from which 
velocity is determined, must contain a sufficient number of cross- 
overs. We use circular patches with radii corresponding to 1° in lon- 
gitude. The patches are located in a regularly spaced grid within the 
ice shelf with intervals of 0.125° in latitude and 0.5° in longitude 
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Fig. 3. Estimation of ice velocity in a test area, (a) Location of test area (red square), (b) Shape of fitting function in a pre-defined velocity space with a maximum deviation of 200 m 
yr _ 1 from a given initial velocity. The dashed and solid arrows show the initial velocity and the estimated velocity. The color lines overlaying the MOA subscene in (c) and (d) show 
surface undulations along the ICESat tracks after relocating with (c) the initial velocity and (d) the estimated velocity. 


over the RIS (3989 assessment patches, roughly 15 km in diameter at 
the latitude of 82°S, -380 crossovers per patch in average). The initial 
velocity fields are linearly interpolated and extrapolated from the 
VELMAP velocity as described in Section 2.3. Fig. 4 shows the new 
ice velocity field over the RIS compared with the initial velocity 
field. The mean value of the deviation from initial velocity is about 
— 48 m yr _ 1 . The largest differences between initial and new velocity 
are observed on the northwestern and southernmost part of the ice 
shelf and the ice front near Roosevelt Island. In the northwestern 
part, our velocity estimation correctly adjusts the initial velocity dis- 
tribution that contains artifacts due to the data gap in VELMAP 
(Fig.l). Negative velocity differences over a large area on the south- 
ernmost part indicate the recent slowdown of the Siple Coast ice 
streams A and B (now Mercer and Whillans ice streams) 
(Bindschadler & Vornberger, 1998; Joughin et al., 2002). The positive 
velocity difference along the ice front near the Roosevelt Island re- 
sults from the velocity discontinuity across a large rift growing after 
the calving of B15 iceberg in March 2000 (Joughin & MacAyeal, 2005). 

4.2. Uncertainty analysis 

The RMS elevation difference of the crossovers at the velocity esti- 
mates v indicates the uncertainty of our scheme, since the ultimate 
goal of the velocity inversion is to minimize the RMS crossover differ- 
ence to the measurement error level. Fig. 5 shows the histogram of 
RMS elevation differences from the 3989 area patches over the RIS. 
The RMS elevation differences of the crossovers at v = v are signifi- 
cantly smaller than those of the crossovers at v = 0 (orbital cross- 
overs) and v = v°, representing the improvement of velocity. 


However, the peak of histogram at v = v (-0.17 m) is greater than 
the peak of measurement error (-0.13 m). We can thus expect that 
the residual error except for the measurement error is -0.11 m ( 

= yl (0.17) 2 — (0.13) 2 ). This residual error results from effects ignored 
in our algorithm, i.e., the elevation changes due to the ice thickness 
change, the velocity variation within a patch (or, equivalently, strain 
of the ice within the patch), and non-linear time variations of 
velocity. 

The formal error of our velocity estimates can be derived from the 
covariance of model parameters, cov(v) = Od(A r A) _1 , where Ais the 
Jacobian at the final solution and Od is the data variance (Menke, 
1989). If the RMS crossover difference of the crossovers at the velocity 
estimates is applied to the data variance Od in each patch, the formal 
errors in average are 4.1 m yr _ 1 and 6.6 m yr _ 1 for the eastward and 
northward components respectively. Fig. 6 shows the distribution of 
formal error. The formal errors of northward velocity component 
are -60% larger than those of eastward ones because the linear fea- 
tures on the ice surface are aligned along the northward direction 
dominantly. Large formal errors along the ice front probably result 
from the rapid basal melts and calving processes. The extremely 
smooth ice surface (e.g. a flat plain) also decreases the signal to 
noise ratio (hence increases the uncertainty of velocity estimates), 
as observed as the pattern with relatively large formal errors aligned 
along the flow direction. Even though the measurements on the 
grounded ice were removed using the grounded ice mask in the 
CATS2008a model, the formal errors are large in some areas near 
the grounding line because of inaccurate grounded ice mask and 
large errors in the tide model. On the other hand, the random noises 
observed in the difference between the initial and estimated 
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Fig. 4. Ice surface velocities derived from ICESat measurements in the Ross Ice Shelf in unit of m yr 1 . The velocity estimates (v x , v y ,v = |v|) are compared to the initial velocities (v Xt 
v° v° = |v°|). 


velocities in Fig. 4 are a little larger than the formal error estimates. 
This means that the actual errors in velocity estimates are slightly 
larger than the formal errors. 

4.3. Validation of velocity estimates 

There are few in situ velocity measurements on the RIS during the 
ICESat campaign period (2003-2009). We compared our results with 
GPS measurements acquired near the ice front in 2005 (Brunt, 2008) 
and InSAR velocity estimates from the Antarctic Mapping Mission 
(AMM-1) in 1997 and the Modified Antarctic Mapping Mission 
(MAMM) in 2000 (Jezek, 2002). Fig. 7a shows the velocity vectors 
at three GPS stations (R13, NASS, and NASC). The velocity vectors 
from the GPS and InSAR measurements are somewhat eastward of 
the initial velocity (from the RIGGS). The directions of the velocity 



Fig. 5. Histograms of RMS crossover elevation differences for intra-campaign cross- 
overs (measurement error), crossovers at zero velocity, crossovers at the initial veloc- 
ity, and crossovers at the estimated velocity. 


vectors obtained in this study are close to those of the GPS and 
InSAR measurements at the NASS and NASC stations. Fig. 7b presents 
the velocity estimates at the three RIGGS stations (J10, I<10, and K1 1 ) 
downstream of Crary Ice Rise, where Joughin et al. (2002) document a 
slowdown of ~69 m yr -1 from 1974 to 1997 based on the velocities 
from RIGGS and InSAR. Our velocity estimates are similar to those 
from InSAR but indicate a slowdown of ~80m yr -1 from 1974 to 
2006. 

5. Discussion 

The velocity estimation from ICESat measurements works well on 
ice shelf surfaces with characteristic topographic features. However, 
under the conditions of extremely smooth ice, ice with rapidly chang- 
ing surface features, or ice with only linear features aligned exactly 
along flow direction, the velocity adjustment we describe would 
fail. If the initial velocity differs too much from the true velocity dur- 
ing the laser measurement period, the minimization of fitting func- 
tion by our non-linear least squares method could converge to an 
erroneous local minimum, resulting in false velocities. Our assump- 
tion that regional ice thickness changes can be ignored fails in areas 
where large and rapid basal melting occurs (e.g. along the ice front). 
The melt rate in the front of RIS is as high as 2.8 ± 1 .0 m yr _ 1 ( surface 
elevation change of 0.31 ±0.1 m yr -1 ) and exponentially decreased 
with the distance from the ice front (Horgan et al., 2011). According 
to their estimates, the surface elevation change due to basal melt is 
able to exceed 0.1 m yr - 1 within the front 20 km. To obtain more ac- 
curate velocity in the area with rapid basal melting, the effect of ice 
thickness change should be removed by an additional parameteriza- 
tion (e.g. linear trend fitting) in the inversion. In addition, significant 
deformation due to strain within an assessment patch would violate 
our uniform velocity assumption within patch, and would lead 
much higher residual errors for the fits. 

The limitations of our method cause a few outliers in the velocity 
estimates as shown in the velocity difference map of Fig. 4. We re- 
move the velocity estimates with the formal errors larger than 13 m 
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Fig. 6. Formal error of the velocity estimated from ICESat measurements: distribution of (a) eastward and (b) northward velocity estimates and their histograms (c, d). 


yr - 1 and then apply 3x3 median filter to generate the final product 
of velocity map (Fig. 8a). The velocity differences at the RIGGS sta- 
tions provide information for the dynamical evolution of RIS during 
last three decades (Fig. 8b). The deceleration of Whillans Ice Stream 
is affecting the velocity structure in the ice shelf far away 
(>300 km) from the grounding lines. In the area enclosed by the 
white line in Fig. 8b, the slowdown is 98 ± 34 m yr -1 (~23 %) and 
the deceleration rate is 3.1 ± 1.1 m yr -2 . These values are reasonable 
though they are slightly lower than earlier estimates of deceleration 
on Whillans Ice Stream (Bindschadler & Vornberger, 1998; Joughin 
et al., 2002). For instance, our estimates of deceleration (~2.5 m 
yr -2 from 1974 to 2006) at J10, K10, and Kll RIGGS stations 
(shown in Section 4.3) are lower than the earlier deceleration rate 
(~3.0 m yr -2 from 1974 to 1997) from Joughin et al. (2002). Howev- 
er, the recent decrease of deceleration rate cannot be assured because 
there exists the uncertainty of deceleration rate (>±0.5 m yr -2 ) due 


to the velocity estimate errors (o- RIGGS = ~15 m yr -1 , o- IC ESat = ~8 m 
yr -1 , and cr InS AR = ~2 m yr -1 at J10, I<10, and Kll RIGGS stations). 

The upcoming polar laser altimeter mission, ICESat-2, will be a 
continuously-operating, photon-counting mission in the 91 -day re- 
peat pattern (Abdalati et al., 2010). The main sensor, the Advanced 
Topographic Laser Altimeter System (ATLAS), will have 6 beams, ar- 
ranged in three pairs, each pair separated by 90 m, with pairs spaced 
3 km apart. With pointing knowledge of around 10 m and repeat- 
track pointing accuracy of a few tens of meters, ICESat-2 will in effect 
provide three swaths of elevation plus absolute slope knowledge, 
with a vertical accuracy of 5 to 10 cm. Applying our technique to 
ICESat-2 data, assuming an excellent tide correction, could yield a de- 
tailed velocity mapping up to 4 times per year, or, after approximately 
2 yrs, a detailed velocity change map with similar accuracy to the 6- 
year map we have derived from ICESat data here. Moreover, the con- 
tinuous operation in the 91 -day repeat pattern, and the three beam- 



Fig. 7. Comparison of the velocity estimates at (a) three GPS stations near the ice front and (b) three RIGGS stations downstream of Crary Ice Rise: black arrows— RIGGS, green ar- 
rows in (a)— GPS, blue arrows— InSAR, and red arrows— this study. The large black circles indicate the patches closest to the stations. The velocities from InSAR were (a) extracted 
from the velocity product of Modified Antarctic Mapping Mission (MAMM) and (b) digitized from Joughin et al. (2002). The background image is the MODIS mosaic of Antarctica. 
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Fig. 8. (a) Ice velocity in the RIS after removing the patches with large formal errors 
(oicESat>13 myr -1 ) and applying the 3x3 median filter and (b) the velocity dif- 
ference between RIGGS measurements and ICESat estimates. The color of dot in 
(b) indicates the velocity difference and the size of dot indicates the combined 

error ^ = ^/of C ESat + °riggs ) contained in the velocity difference. The white line 
encloses the area where large slowdowns are observed. 


pair configuration, will yield ~54 times the number of crossover re- 
gions (i.e. velocity measurement points). A further improvement in 
an ICESat-2-based application of our technique comes from the 
availability of both slope and elevation data continuously along the 
beam-pair swaths. This in effect provides two parameters for the 
RMS fitting routine to optimize, improving accuracy, and leading to 
maps of changing deviatoric strain as well as changing velocity over 
time. 

6. Conclusion 

The ice velocity mapping on ice shelves using the high-accuracy 
and high-resolution laser altimetry can be an alternative scheme to 
investigate the dynamical evolutions of large-scale ice shelves. The 
velocity adjustment minimizing the elevation difference of the cross- 
overs attached on the ice shelf surface achieves accurate velocity esti- 
mates in the presence of surface (topographic) undulations but in the 
absence of time-varying thickness changes or surface deformations. 
Although the inferred errors of velocity estimates are larger than 
those of InSAR application, there are still some advantages. The 
altimetric approach can provide a long-term averaged velocity with 
low costs of computation and data storage whereas the InSAR uses 
relatively short repeats (up to tens of days). Accurate geolocation of 
ICESat data facilitates estimating the absolute motion of ice so that 
the absolute calibration of velocity data, which is a challenging issue 
in the InSAR velocity mapping in Antarctica, is not necessary. 


To overcome the limitation of our method (discussed in 
Section 5), further improvements must be studied, i.e., including the 
vertical velocity (elevation change rate) and regional strain rates 
into the algorithm as additional parameters to be estimated simulta- 
neously. Another potential improvement would be to use patches 
adaptive to glacier boundaries and rifts. It is also recommended that 
observations from more sophisticated satellite altimetry missions 
like namely CryoSat-2 and ICESat-2 are applied to our technique or 
jointly processed. 

The method suggested in this study would work well on the 
Ronne-Filchner Ice Shelf as well as the Ross Ice Shelf but would 
have limited success on small ice shelves or the ice shelves at lower 
latitudes (e.g. Amery, Fimbul, Getz, or Larsen C) because of low 
crossover densities. An extension of our scheme to the simultaneous 
mapping of velocity and vertical elevation change is more important 
for the Ronne-Filchner, Amery or Larsen C ice shelf, because of their 
active basal freeze/melt processes. 
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